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Abstract 

Q-ij We investigate the effect of the bulk contents in the DGP braneworld on the evolution of the 



universe. We find that although the pure DGP model cannot accommodate the transition of the 



effective equation of state of dark energy, once the bulk matter T 5 5 is considered, the modified 
model can realize the w e ff crossing —1. However this transition of the equation of state cannot 
be realized by just considering bulk-brane energy exchange or the GB effect while the bulk matter 
contribution is not included. T| plays the major role in the modified DGP model to have the w e ff 

O ' 

crossing —1 behavior. We show that our model can describe the super-acceleration of our universe 
with the equation of state of the effective dark energy and the Hubble parameter in agreement 

•rH , 

^ ' with observations. 
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I. INTRODUCTION 



The accelerated expansion of our universe is one of the most important discovery in 



Q0 



the last decade [1|, [2j, |3j], having triggered plenty of efforts to understand and explain it. 
This phenomenon is in conflict with our common sense about attractive gravity. Within 
the framework of general relativity, the acceleration is attributed to the mysterious "dark 
energy" existing in our universe. The theoretical nature and origin of this dark energy are a 
source of much debate. Candidates suggested for this dark energy can be classified according 
to the behavior of their respective equation of state w — P/p. The cosmological constant, 
with w = — 1, is located at a centra position among dark energy models both in theoretical 
investigation and in data analysis {4]. In quintessence 5], Chaplygin gas [6] and holographic 
dark energy models ?], w always remains bigger than —1. The phantom models of dark 
energy have w < — 1 8|. Recent more accurate data analysis tells us a dramatic result, 
namely that the time varying dark energy gives a better fit than a cosmological constant 
and in particular, w can cross —1 around z = 0.2 from above to below [9]. Theoretical 
attempts towards the understanding of the w crossing —1 phenomenon have been suggested, 
including the model containing a negative kinetic scalar field and a normal scalar field 10) , 



a single scalar field model llj], interacting holographic dark energy models [12J and others 

An alternative approach which does not need dark energy to explain the late-time ac- 
celeration is motivated by string theory via the brane-world scenarios. In this scenario our 
universe is a 3-d brane embedded in a space-time with extra dimensions. The cosmolog- 
ical evolution on the brane is described by an effective Friedmann equation incorporating 
non-trivially with the effects of the bulk onto the brane. The presence of the 5-d mat- 
ter can interact with the matter contents on the brane and alter the cosmic expansion 
leading to a behavior resembling the dark energy. The cosmic evolution of the Randall- 
Sundrum(RS) _ brane world 14J with energy exchange between brane and bulk has been stud- 



ied 



15 



16, 



In these models, due to the energy exchange between the bulk and 



the brane, the usual energy conservation law on the brane is broken and consequently it was 
found that the equation of state of the effective dark energy can experience the transition 



behavior 



1 i fl 



3- 



In string theory, in addition to the Einstein action, some higher derivative curvature 



2 



terms have been included to derive gravity. The combination of the Einstein-Hilbert and 



Gauss-Bonnet (GB) term constitutes, 
produce second-order field equations 



or 5D spacetimes, the most general Lagrangian to 



20 



21] . The GB correction changes the bulk field 



equations and modifies the braneworld Friedmann equation. It influences the evolution of 
the universe in our brane. Effects of the GB correction on the RS braneworld have been 
studied in 



19 



22j. 



In this paper we are going to concentrate on another braneworld model introduced by 
Dvali, Gabadadze and Porrati (DGP) 23J, where the braneworld is embedded in the flat 
bulk with infinite extra dimensions. Considering that the graviton propagates into the extra 
dimension, and at large scale, gravity can become weaker due to its leakage, the DGP 
model can realize the accelerated expansion naturally. However for the pure DGP model, 
its effective equation of state never goes down to the phantom phase. Our main motivation 
here is to investigate the effects of the bulk contents in the DGP braneworld on the evolution 
of the universe and explore the possibility of the transition of equation of state if there are 
contributions from the bulk-related energy-momentum tensor components which has been 
observed in RS model|15j,lla,ll7|,ll8|,ll9j. The DGP model only with T 5 ° has been investigated 
in 24|. We are going to present a systematic and complete examination of the bulk effects 
including T° and T| terms on DGP model. Besides we will also study the modification on 
the brane evolution due to the GB correction together with bulk related energy-momentum 
tensor components. Influences of the GB correction on the pure DGP braneworld have been 



studied in 125 



261 ]. Although the effects of the GB correction term on the late time universe 



is small, we will see that it still plays an important role in the early time cosmic evolution. 
We will show that in the DGP model the bulk matter contribution T| plays a major role 
in accommodating the transition of equation of state, while the T° and GB correction alone 
cannot present a profile of the w e ff crossing —1 phenomenon found by observations. 

The organization of the paper is the following: in section II we will give out the basic 
equation sets for the DGP model by considering different correction terms respectively. The 
bulk effects due to T| term will be shown in detail in section III. In section IV, we will 
consider the influence of the energy flow T° on the brane universe evolution. Conclusions 
and discussions will be presented in the last section. 
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II. GENERAL EQUATIONS FOR DGP MODEL WITH GB CORRECTION 



The DGP brane model with GB correction starts from the action 

S = / d 5 X^(R - 2A 5 + olLqb) ~ / d 4 x^{R - 2A 4 ) + J d 5 ^g~L EM , (1) 

where k and \i are related to the gravitational constants and the Planck masses for the bulk 
and brane as [27] ]: 

K 2 = 8nG {5) = M 5 " 3 ; /i 2 = 87tG (4) = M 4 " 2 , (2) 

respectively, A5 and A4 are cosmological constants for the bulk and brane. L B m is the 
energy-momentum tensor and Lgb is the GB correction term in the form 

Lgb = R 2 — AR AB R AB + R ABCD R ABCD . (3) 



a is the coefficient of the GB term, which is positive, as required by string theory 20| and 



is general 
model 



y considered to be very small. If we take a = 0, Eq.(l) reproduces the pure DGP 
271 ]. Throughout the paper the capital letter are used to present the 5-d indices, 
while the Greek alphabet is used for 4-d brane case. 
From the action one can obtain the field equation 

Gab + A 5 9ab + 2aH AB = k 2 {T ab - [\{G^ + A 4 ^„) + f^]8{y b )8^8 B }, (4) 

where H AB = RR AB -2R C A R BC - 2R CD R ACBD + R c a de R BC de ~ \qabLgb is the second- 



order Lovelock tensor 2l|, S(yb) comes from the difference between the integration with 4-d 
metric and 5-d metric. 

The energy-momentum tensor on the brane is assumed to be that of a perfect fluid, 

T„ v = (p + pju^Uv + pg^, (5) 

where w M , p and p are the fluid velocity, energy density and pressure, respectively (c = 1 is 
used). The non-zero components related to the fifth dimension in the bulk energy-momentum 
tensor are supposed to be T° and T|, whose role in the accelerated expansion will be studied 
in detail. 

Generally, the metric in 5-d brane cosmology is written as 

ds 2 = -n 2 (t, y)dt 2 + a 2 (t, y)^ ij dxidxj + b 2 (t, y)dy 2 , (6) 



where y stands for the extra dimension orthogonal to the brane, and 7 U is the maximally 
symmetric 3-d tensor. Then yj—g = by/—g, thus 6(yb) = ^jp-. From this metric the Einstein 



b 

equation can be obtained directly. According to Eq.(4) one can obtain the Einstein tensor 
components as [28] 

_ or 2^ ah n 2 a" 
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r 1 a 1 a I; 

$ = 1 . (8) 

n 2 a 2 b 2 a 2 a 2 



The dot denotes a derivative with respect to t, and the prime the derivative with respect to 
y. Without loosing generality, in the brane world scenario, one usually chooses the metric 
function b(t, y) — 1 and n(t, 0) = 1 to simplify the calculation. 

We choose the brane to be located at y = 0, and suppose the metric functions to be 
continuous at this point, but their first derivatives are discontinuous due to the energy- 
momentum tensor distribution on the brane. Furthermore, the geometry is supposed to 
display a ^-symmetry around y = 0, thus a'(0+) = — a'(0_) and n'(0+) = — n'(0_). If one 
integrates the field equation for the tt and ij components at the infinitesimal region near 
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y = 0, only those terms in the metric with a" or n" and the energy-momentum distribution 
on the brane can remain. Then the differences of a' and n' on both sides of the brane, say, 
a'(0 + ) — a'(0_) = 2a'(0 + ) and n'(0 + ) — n'(0_) = 2n'(0 + ) can be obtained. For simplicity, 
throughout the paper we use a' and n' to stand for a'(0 + ) and n'(0 + ), and in the equation 
on the brane all bulk terms are taken with their values at y = + . From G u , Gij, H tt , 
in Eq.(7), we find that a' and n' satisfy the following equations: 

-{-3 + 4a — - 3 - + - } = — /i 2 p - 3 - + - ; 9 
a a z a z a z 2fi z a z a z 

. 2a' r ,,k a' 2 d 2 . 2a' a dn N1 k 2 , 9 k a 2 2an 2d . 

n' + — + Aa[n'{- -—+-+— = —(fi 2 p + — + — + — ). 10 

a a z a z a z a a a 2/j z a z a z a a 

Generally, one can solve these equations and substitute the results of a' and n' into 
the field equation for the ty and yy components to obtain the continuity equation and the 
effective Friedmann equation. But before that, it is helpful to notice that the function $ we 
introduced in Eq.(8) satisfies 

$' = ^(A 5 -T >V-^T 5 Va, (11) 
4 = ^-(A 5 - T|)a 3 a - ~^af, (12) 



as 



where $ = ($ + 2a$ 2 )a 4 . From (12), if T 5 ° = and T| has a proper ansatz, such 
— f flt '> $ can be obtained analytically by an integration with respect to t: 

~ k 2 . 2Fa u+i „ , . 

i = T A ' a -WT7j + c - < 13) 

where C is an integration constant. For the case without GB term, the solution of $ is 
simply 

^ k 2 . 2Fa u C 

$ = — A 5 ; r + —r, (14) 

6 3(4 + i/) a 4 ' V ; 

where the term ^ is usually referred to the dark radiation 29j. For the case with GB 
correction, solutions of $ are, 



ffi ~(4 + .) ± y/(4 + ,) yi + 8a(^ + g ~ jgj) _ 
$ " 4a(4 + z,) ' (15) 

while only one solution with finite a — > limit can be taken. For v > —4 the solution reads 



When a -> 0, Eq.(16) goes back to Eq.(14). 

From the definition of <£> in Eq.(8), we see that a' can be expressed in terms of a once $ 
is obtained. Integrating the equation of tt component around y = and substituting a' in 
terms of <J> and d, we can finally arrive at the equation for the Hubble parameter H(t) = ^, 

,^ + 1 _ + + * + | )P = + * _ + Al)P , (17) 

where r = h 2 j p 2 is the DGP crossover radius. 

For the sake of simplicity and clarity in the following discussion, we give the simplifications 
we are going to use. We will neglect the cosmological constant, A 5 = A 4 = 0, since the effect 
of the cosmological constant can be included in p and p. We will apply our discussion to the 
flat universe with k — 0. Besides, we will employ dimensionless notations in the following 
calculation by defining 

H 2 



x = 



u 



y = 



n = 



ao _ 1 
a 

p 2 p 



H 2 r 2 ' 

8H 2 
m = — -a, 



~ H 2 {A + u) ' 

X - a o p 
~ H 2 ' 



ZM a 

where a and if are the present values of the scale factor and the Hubble parameter, z is 

2 

the redshift and u/x = fj^ is the proportion of matter in the total effective energy density. 
Using dimensionless notations, the expression for the solution of $ becomes 

* = 2H„- X{l + ^ + M{l + ^, (19) 

for the case without the GB correction; if the GB term is included, it reads 

-1 + Jl + 2m(-X(l + z)~ v + Mil + zY) 

$ = 2H y ; — . (20) 

3m 



The equation for H 2 , Eq.(17), turns into an equation for x [26| 

An(x - 2y)[l + m(x + y)] 2 = (x - u) 2 . (21) 

If T° is nonzero, $ cannot be solved analytically and we can not take advantage of 
the simplicity discussed above. To obtain the equation for H 2 , we need to substitute the 
solutions of a' and n' into the equation of yy component, and obtain H(t) through onerous 
calculations. The acceleration of the scale factor can be written as a = a(H 2 + H). By 
choosing the proper ansatz of T5 and expressing the result of p(t) as a function of a(t), 



we can generally obtain the equation as a nonlinear ordinary differential equation of H 



m, 



(t) 



combining with the unknown function a{t). To solve such a problem in RS model 
the authors introduced new effective fields related to H{t) 2 and separated the equation into 
two equations, both of which are solvable separately. But for the DGP model, due to the 
extra 4-d intrinsic curvature terms, the highest order of H(t) is 4, rather than 2 in the RS 
model. With the GB correction, the order goes up to 8. In Ref. 24j when just nonzero T° 
was included in the pure DGP model, the author solved the problem by introducing the 
concept of "fix point" and setting p(t) and the auxiliary field time-independent. Generally, 
we do not hope to obtain any analytical solution for such a nonlinear ordinary differential 
equation. We will count more on the numerical calculations. Considering nonzero T 5 ° and 
T| components, our problem is general and complicated. We will present a general way to 
solve the problem. 

We have two time- dependent functions, H(t) and a(t), in the same equation. Consider- 
ing that in the big-bang cosmology the flat universe is expanding monotonically, a(t) is a 
monotonic function of t, and H(t) can be written as H(a). For the convenience we will use 
the dimensionless redshift z and write H(t) as 

H{t) = 2 dz~- (22) 

Substituting all the dimensionless notations into the equation of the yy component and 
expressing the results until the linear order in a, the equation for H(t), or equally, x(z), is 

= — mi 4 + m[2An + u + (1 + z)x']x 3 

+ [— 16n + 4877m 2 — Ylmnu + 3mu 2 — 3m(6n + u)(l + z)x'\x 2 

+ {64?7 2 + 87777 — 12mnu 2 — 5m« 3 + [877(1 — 3777,77) + 3mu(8n + 77)] (1 + z)x'}x 



8 



+ [8nu 2 + 2mu 4 - {16n 2 + u[8n + mu(Qn + «)]}(1 + z)x'] 

+ 32 " 2(1 3 +Z) ""x, (23) 

where the prime here is the derivative with respect to z, x and u are functions of z, and we 
have taken T| = ^a v ' . It is to be noted that since we don't need to analytically integrate 
Tf with the term a 3 d as did in Eq.(12), we can in principle use any form of T| as a function 
of a. Eq.(23) is a nonlinear differential equation of x(z). 

From the equation of the ty component, assuming b(t,y) = 1, we get 

3(l + 4a$)(----) = T 05 . (24) 
an a 

Taking the value of each term in this equation at y = 0+ and substituting the solution of a' 
and n', we can find that the left hand side of this equation is simply \{p + 3if (p +£>)). If 
T° = 0, it is the conservation of energy on the brane. If T° ^ 0, it acts as the energy flow 
between the brane and the bulk, 

p + 3H(p + p) = -2T°. (25) 

Here T 5 ° has a sign difference as compared to T 05 due to the metric term gtt\ y =o = —n(t, 0) = 
-1. 

Setting the ansatz, T 5 ° = fHa s , the equation for p can be solved analytically. Expressing 
pit) as p(a), we have 

jg + 3 ( 1 + ^ + 2/q'- 1 = 0, (26) 
da a 

with the solution 

P = a-^C, - 2/a ' , (27) 
3 + 3w + s 

where C\ is an integration constant. For the cold matter on the brane w — 0, the first term 
on the right- hand- side is proportional to a~ 3 and the second term could be attributed to the 
effective dark energy. Eq.(27) can be expressed by using the dimensionless notation u(z), 

u(z) = P(l + z)- s + n m0 (l + z) 3 , (28) 

where P = ~ 3(3+3°^ /; an d ^ m o — sHfia? = wfP° ^ s ^ ne P resen t ratio of conservative matter 
in the total energy density of the universe, p is the density of the conservative matter today. 
There is a strong constraint on the value of dimensionless parameter P. Since the matter 
portion of the total energy density should be in the range [0,1], thus < ^ < 1. Here we 



will use the fitting results on the WMAP data [30J] and take f2 m0 = 0.28 in our calculation. 
Of course, that fitting is from a different model, but the generally accepted values of f2 m0 are 
all close to this value, and the small variation of this value will not change the qualitative 
conclusion of our calculation. At the present moment z = 0, we have —0.28 < P < 0.72. 
Any solution with P out of this range is physically unreasonable. 



To describe the effect of the effective dark energy, we can define the effective equation of 



state 



3l|: 



, 1 d\n5H 2 

w( ^" s - 1 + (29) 

where 5H 2 = H(z) 2 — f2 m0 (l + z ) 3 Hq. w e ff can be expressed by using the dimensionless 
parameters as 

W{ZUS 3x(z) - 3JW1 + zf ■ (30) 

The subscript eff indicates that the effect similar to the dark energy on the brane comes 
from the bulk contribution. Another important quantity is the deceleration parameter q 



da 1 3 , s ^mo(l + z 

which will also be used in the following discussion of the expansion of our universe. 



q = = 2 + 2 W{z)e " {1 ~ x(z) } ' (31) 



III. CALCULATION AND DISCUSSION WITHOUT CONSIDERING THE T 5 ° 
TERM 

For the case without T 5 ° term, we can solve the equation Eq.(21) by substituting $ from 
Eq.(20) and u from Eq.(28) but with P = 0. We can put the result into Eq.(30) and Eq.(31) 
to examine the behavior of the effective dark energy. 

In Eq.(21), if there is no GB correction, the equation for x is quadratic, thus we have two 
solutions. For the case $ = 0, they recover the two branches of pure DGP model, DGP(+) 



and DGP(-): x = 2n + (1 + z) 3 Q m0 ± 2jn 2 + n(l + zfVL m0 . Since only the DGP(+) 



251 ] . we will concentrate our discussion on 



'mO 

solution has late-time self-accelerating behavior 
this solution. When GB correction is added, Eq.(21) becomes a cubic equation and has 
three roots for x, two of which correspond to DGP(+) and DGP(-) in a — > limit, and the 
third one diverges when a — > 0. For comparison, we also study the solution with DGP(+) 
limitation for that case in this work. 



10 



To compare our model description on the evolution of the universe with the observation, 
we have several constraints to meet. At present we have x(z = 0) = H(t = 0) 2 /Hq = 1. 
7 or the effective equation of state, we require w e jf(z = 0.2) = —1 and w e ff(z = 0) = —1.06 
9J. We will use these constraints to refine our model parameters. Assuming If = 0, we 
have parameters such as n (corresponding to the DGP crossover radius), m (corresponding 
to the GB correction), M (corresponding to the dark radiation), X and v (relating to T| 
form). We will focus on whether we can accommodate the w e ff crossing —1 by including 
bulk related energy-momentum tensor and the GB correction, which cannot be realized in 
the pure DGP model. 

Since we have the parameter v in the exponential, the equations are not polynomial. We 
will use FindRoot in our calculation which will raise the problem about the choice of the 
initial values. To avoid the possibility of failing to find the solution, we will try different 
initial values in the FindRoot. We find that the solution is not strongly dependent on the 
initial values so that we are confident that our results are almost the whole collection of all 
possible solutions to the equations we deal with. 

To see the consistency of our results with observation, we will plot the Hubble parameter 
and compare with the observational H(z) data as shown in Table 1. It is interesting to note 
that all the cases which can accommodate the equation of state transition can fit well the 
observational H(z) data. Remembering that the u> e // crossing —1 was observed in the SNIa 
data fitting containing an integration effect in the luminosity distance, while the Hubble 
parameter does not suffer from this integrated over effect, the Hubble parameter data can 
present a complementary and consistent check for our model. 

Now we list out our results step by step. In the following results, all the numbers obtained 
in numerical calculation are expressed only with 3 digits after decimal unless for the cases 
where more digits are necessary to be shown. 

1. DGP+M 



z 


0.09 


0.17 


0.27 


0.40 


0.88 


1.30 


1.43 


1.53 


1.75 


H(z) (km/s/Mpc) 


69 


83 


70 


87 


117 


168 


177 


140 


202 


la uncertainty 


±12 


±8.3 


±14 


±17.4 


±23.4 


±13.4 


±14.2 


±14 


±40.4 



TABLE I: The observational H (z) data 



32 



33] 
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In this step we consider the DGP model with the dark radiation, where M denotes the 
dark radiation. Now we have two free parameters n and M and we will use two constraints: 
x(0) = 1 and w e ff (0.2) = — 1 to see whether the dark radiation can help to realize the w e ff 
crossing — 1. Actually n and M can be solved by using these constraints as n = 0.157 and 
M = 0.263. But using these values of n and M, the w e ff behavior is not good. w e ff crosses 
— 1 at z = 0.2 from below to up as shown in Fig.l, and the present value is w e ff = —0.950. 
This is not in consistent with the observation, especially the transition behavior. 




FIG. 1: The w e ff curve as function of z in the case DGP+M. The behavior of w e ff is bad since it 
crossed — 1 at z = 0.2 from below to above. 

2. DGP+GB+M 

Including the GB correction, we have three free parameters now, such as n, m and M. 
If we apply all three constraints (x(0) = 1, iu e //(0.2) = —1 and w e ff(0) = —1.06), the 
solution is complex (n = -0.030 + 0.093i, m = -0.002 - 0.007i and M = 6.340 + 6.244i). 
If we apply only two constraints (x(0) = 1 and w e ff(0.2) = —1) and search n in a big 
range 0.001 < n < 5, the similar result to that in case 1 appears: w e ff crosses —1 from 
below to above and m is negative. One solution is shown in Fig. 2a. We note that there is a 
singularity about z = 1.246 in the curve of w e ff. This singularity comes from the definition 

12 



in Eq.(29), we see that when 5H 2 = H(z) 2 — f2 m o(l + z) 3 Hq < 0, the logarithm is not well 
defined, but one can still calculate the w e ff through the simplified expression in Eq.(30). In 
Fig. 2b we show the relation between H 2 /Hq and the matter component f2 m0 (l + .z) 3 . We see 
that beyond the redshift z = 1.246 the matter component is overweight, so 5H 2 < 0, which 
means the effective dark energy component is negative. This is obviously unreasonable, and 
at least it shows that the model fails in explaining the universe before that redshift. In this 
work we concentrate on those solutions with w e ff(z) free of singularity. 




FIG. 2: In the case DGP+GB+M, w e ff and q as functions of z are shown in Fig. 2a. The behavior 
is not favored as in Fig.l. Singularity is observed at z = 1.246. In Fig. 2b, relation between H 2 /Hq 
and the matter component O m0 (l + ^) 3 is shown. When z > 1.246, 5H 2 = H(z) 2 — fi m o(l + z) 3 H^ 
becomes negative, which breaks the definition of w e ff i n Eq.(29). 

3. DGP+Tf 

Now we include the bulk related energy-momentum tensor T|. In this case we have three 
free parameters (n, X and v) and we are going to employ all three constraints (x(0) = 1, 
w eff (0.2) = -1 and w eff (0) = -1.06). We can find the solution n = 0.046, X = 2.729 and 
v = 0.948. The curves of w e ff, q and H v.s. redshift z are shown in Fig. 3a and Fig. 3b 



respectively. In plotting Fig.3b, we have used H = 72km/s/Mpc 34j. It is interesting to 
find that parameters adjusted to meet the requirement of w e ff crossing —1 and its value at 
the present moment automatically fit well to the H(z) data. Recalling that the transition 



13 



behavior of w e ff results from the SN data analysis containing integration in the luminosity 
distance, while the Hubble parameter is not integrated over, which persists fine structure 
highly degenerated in the luminosity distance, the simultaneous satisfaction of the vo e ff 
behavior and the H(z) data gives complementary and consistent check of the viability of 
our model. 



0.5 




Fig.3a Fig.3b 



FIG. 3: w eff and q v.s. z (Fig.3a) and H(z) curve (Fig.3b) in DGP+T 5 5 case. We see that H(z) 
curve fits the data quite well. 

4. DGP+GB+T| 

Here we include the GB correction based on the case discussed above. We now have 
four free parameters n, m, X and v. Employing constraints (x(0) = 1, u> e //(0.2) = —1 and 
w e ff(0) = —1.06) and searching through 0.001 < n < 5, we see that m, X and v can either 
be negative or positive, but the latter two are always with the same sign, m decreases with 
the increase of n and can be positive only when n < 0.05. The w e ff curve has no singularity 
for positive m, but contains singularity when m < (few solutions without singularity have 
been found with negative m, but then u becomes smaller than —4, which conflicts with our 
simplification assumption discussed above). We are interested in the positive m, since GB 
correction coefficient a is always positive, the singularity- free curves of w e ff(z) and q(z) with 
positive solution m = 0.025 are shown in Fig.4a. When the GB correction is considered, 
w e ff appears more different at larger z if compared to the result without the GB correction. 
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Since the GB effect is only important in the early universe, its stronger modification to the 
w e ff at bigger redshift is natural. In Fig.4b, we plotted the H(z) curve by using the same 
adjusted parameter from the constraints on the equation of state, and we see again that in 
the model when the w e ff requirement is met, the H parameter automatically fits the data, 
which gives the consistent check of the model. 




FIG. 4: w e ff(z) and q(z) curves (Fig.4a) and H(z) curve (Fig.4b) in DGP+GB+T| case. Compar- 
isons with the results without GB corrections have been shown. Differences from the case without 
the GB correction become bigger at higher redshift. 

5. DGP+T|+M 

Based on case 3, we include the dark radiation contribution. We now have four param- 
eters, namely n, X, v and M. Employing three constraints (x(0) = 1, iu e //(0.2) = —1 and 
w e ff(0) = —1.06), and searching n in the range 0.001 < n < 5, we find that the solutions 
exist only when n < 0.36. We see that M can either be positive or negative: for the negative 
M, w e ff never decreases with the increase of z within z < 5; while for the positive M, w e ff 
drops at large z, and singularities of w e ff(z) and q(z) appear within z < 5. Pictures of these 
two cases are shown in Fig.5a and Fig.5b. It is also observed that with the increase of n, 
the positive value of M becomes bigger and the singularity appears at smaller z. 

6. DGP+GB+Tf+M 

Now we have all five parameters: n, m, X, v and M and three constraints to be used 
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FIG. 5: We show in Fig.5a the w e jf(z) and q(z) curves in DGP+T| +M case, where M is negative 
and the curves are singularity free; in Fig.5b, M is positive and the curves have a singularity at 
z = 0.936. 

(x(0) = 1, w e ff(0.2) = — 1 and w e ff(0) = —1.06). First, we try to set m small, e.g. 1CT 9 , 
and choose one solution obtained in case 5, but the FindRoot command cannot help to get 
a solution to meet all three constraints, even though the initial parameters are set closely 
to those in case 5. This means that the solution is sensitive to a, though a (or m) is small. 
Next we search the nonsingular solutions within the parameters' ranges 0.001 < m < 0.02, 
0.01 < n < 0.1, the curves are not quite different from what we have seen in case 3 and case 
5. We show one of the solutions in Fig. 6a and Fig.6b. 



IV. CALCULATION AND ANALYSIS CONSIDERING THE T 5 ° TERM 

The study with nonzero T 5 ° becomes more complicated. In order to solve the effective 
Friedmann equation, or equivalently, to obtain x(z), we have to use the differential equation 
(23) with the full expression of u in Eq.(28). The free parameters we have now are n, m, 
X, z/, P and s. Here we do not have the explicit dark radiation term M, since we do not 
make the integration to get $. In the calculation, we have to solve the differential equation, 
where the boundary condition x(z = 0) = 1 is needed. The effect of the dark radiation is 
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FIG. 6: w ef f(z) and q(z) curves (Fig.6a) and H(z) curve (Fig.6b) in DGP+GB+T|+M case. 

reflected in the boundary condition of x. Solving equations when T° = 0, the result is the 
same as the case when dark radiation term M appeared in the last section. We will still use 
two constraints on w e ff, which are w e //(0.2) = —1 and u> e //(0) = —1.06. 

Eq.(23) is a differential equation of x(z) ) where free parameters are involved. More ef- 
forts are needed to solve the equation numerically. Here we try to employ the self-consistent 
method. First we get x(z) solved with the chosen initial values of parameters, then we sub- 
stitute the numerical result into the expression of w e ff(z). The term in the expression 
of w e ff can be replaced by the function of x(z) using Eq.(23), thus we can write w e ff into 
a function of x(z). Substituting the solution of x(z), w e ff becomes a function of z and we 
can solve the parameters with the constraints w e //(0.2) = —1 and u> e //(0) = —1.06. If the 
result is not consistent, we substitute the results back to x(z) as new initial values until 
convergence is finally arrived. A proper choice of the initial parameters is crucial to obtain 
a convergent result, we usually do the iteration with many different choices within quite a 
large reasonable parameter space. 

Now we show the results we have obtained. First, to show the consistent and efficient of 
our numerical calculation, we turn off the contribution of T° to recover corresponding cases 
in section III. 

7. DGP+T| (employing Eq.(23) to solve the problem numerically) 
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We have three parameters n, X and v in this case. Searching within ranges — 100 < X < 
100 and -4 < v < 100 by setting n = 0.01, we can find the solution (X = 79.824, v = 0.719), 
which is singularity free at least for z < 5. This result can be compared with that in case 5 
(DGP+T|+M). Here the dimensionless parameter for T| term is X, which differs from X 
employed in case 5 with a factor 4 + v as shown in Eq.(18). Taking this into account, the 
solution for X is 16.917, which is just the result in case 5 with n = 0.01, where the solution 
in case 5 reads X = 16.917, v = 0.719 and with the additional parameter M = —1.023, 
see Fig. 5a. They coincide, although they are obtained by completely different methods. 
This shows the correctness of the self-consistent method we used and also demonstrates the 
equivalence of the boundary condition in differential equation and the extra freedom of the 
integration constant. 

8. DGP+GB+T| (employing Eq.(23) to solve the problem numerically) 

We have now four parameters n, m, X and v. Within ranges — 100 < X < 100, < z/ < 5 
by setting n = 0.05 and m = 0.007 (which are the values used in case 6 in Fig. 6a and Fig.6b), 
we can find the solution X = 11.93120 and v = 0.98898. Considering the difference between 
the dimensionless notations, this corresponds to X = 2.39151 and v = 0.98898, which is a 
bit different from those directly obtained in case 6 as X = 2.39146, v = 0.98906. This small 
difference is due to the approximation we have taken in Eq.(23), where the expansion on a 
is kept only to the linear order. So when q ^ 0, the equation systems in case 6 and case 8 
are not exactly the same. We can see that the difference between the curves in Fig. 7a, 7b 
and those in Fig. 6a, 6b lies in large z region, which shows that the effect of GB correction 
is important in the large redshift era. 

To demonstrate more explicitly the effect of GB term, we calculate in this case with 
different fixed values of m. We shut down the freedom of v by setting v — 1 in order to show 
the effect of GB term more clearly. The results are shown in Fig. 8, where we see clearly that 
the GB term only changes the property in the early universe. 

9. DGP+T 5 ° 

We turn on the T° effect to consider the energy exchange between the bulk and the brane. 
In this case free parameters are n, P and s. Searching in the range 0.01 < n < 0.2 with 
-0.28 < P < 0.72 and -100 < s < 100 as initial tries, we find that for n < 0.15, solutions 
can be found, but all require P < —0.28, which is forbidden as we discussed. This tells us 
that T 5 ° alone with the simple form T 5 ° = fHa s cannot lead to the expected behavior of the 
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FIG. 7: w e ff(z) and q(z) curves (Fig.7a) and H(z) curve (Fig.7b) in DGP+GB+T 5 5 case. This 
result is calculated from the equation considering the Tg term, which is not quite different from 
Fig. 6 of case 6, except the curves lie a little higher than the curves in Fig. 6a at large z. 

effective equation of state. 

10. DGP+GB+T 5 ° 

We have four parameters n, m, P and s, one more than those in the previous case. But 
after searching in ranges 0.001 < n < 0.1, 0.01 < m < 0.9, with -0.28 < P < 0.72 and 
— 100 < s < 100 as initial values, all solutions which can be found needs P < —0.28. The 
value of P increases with the increase of m and the decrease of n, but it can only go up to 
—0.425 when n = 0.001 and m = 0.9, which is almost the most favored parameter set within 
acceptable ranges for m and n. Actually m is related to the GB correction which should 
be small. Thus introducing one more free parameter, the GB correction, cannot change the 
unfavored result in case 9. 

11. DGP+T|+T 5 ° 

We contain now five parameters in total, e.g., n, X, u, P and s, but for simplicity, we 
keep only three parameters (n, P and X) free, by setting v — 0.7, s — 1. The reason of 
choosing v = 0.7 is because in case 7 the solutions give the value of v around 0.5 ~ 1. We 
find that the value of n cannot be too big, otherwise the curve will have singularity at very 
small z. For small n, we can find solutions without singularity, e.g., n = 0.01, P = 0.020 and 
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1 2 z 3 4 5 

FIG. 8: w e ff(z) curve in DGP+GB+Tf case. These curves correspond to m = 0.01, m = 0.001, 
m = 0.0001 and m = 0.00001 respectively. 

X = 78.535. We show the proportion of different components as functions of z in Fig. 9, in 
which ^t m atter + ^darkenergy = 1, ^matter is the total matter as obtained from the differential 
equation (26); while ^effective is the effective dark energy proportion including the energy 
exchange effect between the bulk and brane. Since the energy exchange effect is very small 
(P is small), the difference between VL e ff ective and Qdarkenergy can basically be neglected. The 
behavior of w e ff{z), q(z) and H(z) are shown in Fig. 10a and Fig. 10b. Different from the 
case 9, we see that when the bulk matter T| is considered, the modified DGP model allows 
the w e ff crossing —1 and is consistent with H(z) data. 
12. DGP+GB+T|+T 5 ° 

This is the most general case in our discussion, where we have all parameters of our 
model: n, m, X, v, P and s. To simplify the calculation, we fix v = 0.5 and s — 1. 
Searching in parameters' ranges 0.01 < n < 10, 0.001 < m < 1, -0.28 < P < 0.72 and 
— 1000 < X < 1000, we found that when n becomes larger, the solution becomes worse, 
either the curves have very bad shapes or the value of P lies far away from the acceptable 
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1 2 z 3 4 5 

FIG. 9: Different components as functions of z in DGP +T| +Tf case. The solid line and the long 
dashed line, representing the "real" matter component and dark energy component respectively. 
The short dashed line shows the remainder of the total energy density after subtracting the con- 
served matter fi m o(l + z) 3 Hq / H (z) 2 , which acts as the effective dark energy where the energy 
exchange was considered. Since P is quite small, the effect of energy exchange is negligible, curves 
of ^effective and ^darkenergy be almost together. 

range. It is found that generally n should not be larger than 0.1. For small n, we can find 
the solution such as n = 0.001, m = 0.01, P = 0.166, X = 895.044, whose corresponding 
curves of w e ff(z), q(z) and H(z) are shown in Fig. 11a and Fig.llb respectively. 

V. CONCLUSIONS AND DISCUSSIONS 

In this work we have generalized the DGP braneworld by including bulk matter content, 
bulk-brane energy exchange and adding the GB curvature correction term in the bulk action. 
We have investigated the effects of the bulk contents and the GB correction on the evolution 
of the universe. We have found that although the pure DGP model cannot accommodate the 
transition of the equation of state as indicated by recent observation, once the bulk matter 
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FIG. 10: w ef f(z) and q(z) curves (Fig. 10a) and H{z) curve (Fig. 10b) in DGP+T 5 5 +T 5 ° case. 




FIG. 11: w eff (z) and q(z) curves (Fig.lla) and H(z) curve (Fig.llb) in DGP+GB+T 5 5 +T 5 ° case. 

T| is considered, the modified model can accommodate the w e ff crossing —1. However this 
transition of the equation of state cannot be realized by just considering bulk-brane energy 
exchange or the GB effect but without the bulk matter contribution. Thus T| plays the 
major role in the modified DGP model to have the w e ff crossing —1 behavior. The GB 
term can have little influence on the late time behavior of the universe, it gives modification 
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to the equation of state at big redshift. This is because of the fact that the GB correction 
arises from the high energy theory, being negligible in our present cold universe. Besides 
the w e ff crossing behavior, our model can describe the Hubble parameter consistently with 
observation. 

In our parameter space there is a generally favored range n < 0.1, which is crucial to have 
singularity free behavior in the equation of state. From the definition n = this range 
of n requires that the crossover factor obeys r > 3.16HQ 1 , which is bigger than the lower 
bound just due to the GB correction [25l. l35j| . Since P and m are related, proper P requires a 
bit bigger value of m. The permitted value P is small. Its sign corresponds to the direction 
of the energy flow. The results we show previously in case 11 and case 12 have positive P 
standing for the influx of energy, which is considered reasonable as the explanation of the 
accelerating expansion of the universe for the cosmology without extra dimension. But in 
the brane cosmology, since we have shown with our results that the T| term dominates the 
effective equation of state behavior, there is no big difference whether the energy flows into 
or out of the brane, and in fact the solutions with negative P have also been obtained, which 
have similar behavior to those shown here. 

From our result we see that the T° term plays little effect in the transition of equation of 
state, this could be due to the choice of the ansatz. A more general form of T® can make the 
calculation more difficult, since the numerical solution rather than the analytical solution 
of p from Eq. (25) will bring more difficulties in the following calculations. The solution of 
p(a) should be substituted into the function of H(z) after it is expressed in dimensionless 
notation. But in principle this is not impossible. Our method to numerically solve the 
nonlinear differential equation of H(z) 2 supplies a general way to deal with such problem, 
and it can relax the assumption form for T° and T|. We expect to see the influence of a 
more general form of T° on the behavior of the equation of state of effective dark energy. 
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